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Abstract 

We study the existence and stabihty of phaselocked patterns and ampUtude death states in a 
closed chain of delay coupled identical limit cycle oscillators that are near a supercritical Hopf 
bifurcation. The coupling is limited to nearest neighbors and is linear. We analyze a model set of 
discrete dynamical equations using the method of plane waves. The resultant dispersion relation, 
which is valid for any arbitrary number of oscillators, displays important differences from similar 
relations obtained from continuum models. We discuss the general characteristics of the equilibrium 
states including their dependencies on various system parameters. We next carry out a detailed 
linear stability investigation of these states in order to delineate their actual existence regions and 
to determine their parametric dependence on time delay. Time delay is found to expand the range 
of possible phaselocked patterns and to contribute favorably toward their stability. The amplitude 
death state is studied in the parameter space of time delay and coupling strength. It is shown that 
death island regions can exist for any number of oscillators in the presence of finite time delay. 
A particularly interesting result is that the size of an island is independent of N when A'^ is even 
but is a decreasing function of when N is odd. 
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I. INTRODUCTION 



The emergence of self-organized patterns is a common feature of nonequilibrium systems 
undergoing phase transitions. Examples of such patterns abound in nature ranging from in- 
tricate mosaic design on a butterfly wing, vortex swirls in a turbulent stream to synchronous 
flashing of fireflies [ij 



The study of pattern formation and self-synchronization using simple 
mathematical models consisting of coupled differential equations has therefore been an active 
area of research spanning many scientific disciplines. For example such models have been 
employed to examine the collective output of arrays of lasers j2[0], coupled magnetrons 0, 
Josephson-j unctions j^, coupled chemical reactors [1, 0| and electronic circuits Rings of 
oscillators have also been useful for modeling biological oscillations |^ and simulation 
of phase relations between various animal gaits These simple discrete models (usually 
composed of coupled limit cycle oscillators) are amenable to some straightforward analysis 
and also direct numerical solutions particularly when the number of oscillators is small. 
When the number of oscillators becomes very large it is often convenient to go over to the 
continuum limit of infinite oscillators and model the system with a partial differential equa- 
tion. Examples of the latter are reaction-diffusion type mathematical equations, such as the 
time dependent complex Ginzburg-Landau equation (CGLE), that have been widely used 
to simulate the dynamics of pattern formation in fluids and other continuum systems jl^ . 
One of the most popular of the discrete models is the 'phase only model' that arose from the 
pioneering contributions of Winfree |i9|] and Kuramoto [13., .14]. This model is valid as long as 
the coupling between the oscillators is assumed to be weak so that amplitude variations can 
be neglected. When the weak coupling approximation is relaxed and amplitude variations 
are retained the system is found to admit new collective states such as the amplitude death 
state 



- This strong coupling model 
lllSHil. The importance of 



chaos and multirhythmicity [U, 
has received a great deal of attention in recent times [Ij 
time delay in coupled oscillator systems and its effect on their collective dynamics has been 
recognized for a long time. In particular there have been a large number of studies devoted 
to the study of time delay coupled 'phase only' models 111 111, Q H 111, H H H H llj . 
More recently the strong coupling model has also been generalized to include propagation 
time delays H, 11 37, 3^. A number of interesting results regarding collective oscillations 
have emerged from t his g eneralized time delayed model some of which have been experimen- 
tally verified IH H 

Our present work is devoted to further investigations and understanding of the generalized 
time delay coupled model. Most of the past investigations on this model have been restricted 
to collective states emerging from global coupling among the oscillators, also known as the 
mean field approximation. We have chosen to study here the equilibrium and stability 
of collective modes emerging from a time delay coupled system of identical oscillators in 
which the coupling is local and restricted to nearest neighbors. The choice of this model is 
motivated by several considerations. From a basic studies point of view the model offers 
us an opportunity to compare and contrast the properties of the collective states of the 
local coupling geometry to that of the globally coupled states. The model further permits 
a parametric study of the influence of the number of oscillators N on the stability and 
existence properties of the equilibrium states. In the limit of large the model also has 
close resemblance to the continuum CGLE and it is of interest to delineate the important 
differences in the collective properties of the two systems. From a practical point of view our 
model has great relevance for many physical systems including coupled magnetrons, study of 
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collective phenomena in excitable media, coupled laser systems, and small-world networks. 
For example, an interesting study on the role of the geometry of the coupled magnetrons on 
their rate of phaselocking was made by Benford et al. Q who experimentally demonstrated 
production of higher microwave power at Giga Watt levels through phaselocking. Our 
analysis of nearest neighbor coupled oscillator array could provide important clues on not 
only the effect of circular geometry but also on the effects of time delay in such a geometry. 
A ring of coupled oscillators is also commonly encountered in a variety of biological clocks 
P]. For simplicity of analysis most mathematical studies of these rings have been confined 
to an examination of their phase evolutions. Our model analysis can provide an extended 
understanding of these biological clocks with the inclusion of time delays and amplitude 
effects. Another interesting and more recent potential application of our work is in the 
rapidly developing area of small-world networks where the simplest studied configuration is 
often a ring of oscillators with nearest neighbor connections dominating i^l- A recent study 
has explored the role of time delays in a ring of connected oscillators with a special focus on 
small- world networks js^]- We believe that our present study of the ring configuration could 
help provide further insights in this direction with the incorporation of amplitude effects. 

The existence of equilibrium states of locally coupled oscillator systems has been studied 
in the past using group theoretic methods [i^ 3|. We use a plane wave method and 
investigate both the existence and stability regions of the equilibrium states of the coupled 
system. We derive a general dispersion relation that can be used to determine the equilibrium 
states of any number of coupled oscillators. The method is further extended to a linear 
stabihty analysis of these equilibrium states. Using this technique we dehneate the existence 
regions of the various equilibrium phase-locked patterns of the coupled system. In the limit 
of small and in the absence of time delay our results agree with past findings, but we 
also find new equilibrium states that have not been noticed before. With time delay the 
existence regions change in a significant manner and in some of the regimes we observe 
multirhythmicity. Finally we also examine the nature of the death state in our system in 
the presence of time delay. Although time delay induced death state in locall y c oupled 
identical oscillators has been numerically observed in one of our past studies [35| there 
has been no systematic analytic investigation of this phenomenon. We carry out such an 
analysis here and delineate the existence of the death region as a function of the time delay 
and coupling strength parameters. Our analysis yields an interesting result in that the size 
of the death island is found to be independent of A^ (the number of oscillators) when A^ is 
even but is a decreasing function of A^ when A^ is odd. 

The paper is organized as follows. In the next section we present the model equations and 
briefly discuss their similarities and differences with the continuum model CGL equation. 
Using the plane wave method we next derive the dispersion relation in Section IIIII and 
proceed to use it to delineate the possible existence regions of the phase-locked patterns 
both in the presence and absence of time delay. Section IIVI is devoted to a linear stability 
analysis of these states in order to identify their actual existence domain. In Section IVI we 
briefly look at the stability of the origin, that is the existence domain of the death state, as a 
function of the time delay parameter and the coupling strength. Our results are summarized 
and their implications discussed in the concluding Section IVII 
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II. MODEL EQUATIONS 



We consider a one dimensional closed chain of delay coupled identical limit cycle oscilla- 
tors that are close to a supercritical Hopf bifurcation. Assuming the coupling to be linear 
and of the nearest neighbor kind, we can write down the following set of model equations 
to describe the time evolution of the oscillator states, 

^ = (1 + tcuo- I + K[i;^+iit - r) - + ir[V^,-i(t - r) - (1) 

j = 1, . . . ,N, where ipj{= Xj + iYj) is the complex amplitude of the j^^ oscillator, K > 
is the coupling strength, r > is a fixed time delay, and c^o is the natural frequency of 
the uncoupled oscillators. The basic nonlinear oscillator we have chosen here is simply the 
normal form of the supercritical Hopf bifurcation equation truncated to the third order. It 
is also known as the Stuart-Landau oscillator and has been extensively used in the past as a 
model equation for studying nonlinear phenomena in fluids, lasers, and Josephson junctions. 
For example, in the absence of time delay (r = 0), the above set of equations ((T)) have 
been used to study the coupled wakes arising from a flow behind an array of equally spaced 
parallel cylinders |4^. Another interesting connection to past work is the continuum limit 
for r = 0, where one can introduce a lattice spacing a, set ipj = V'(io) in ©; rescale K as 
K/a^ and then let a — >■ with ja x. This reduces Eqs. for r = to the well known 
complex Ginzburg Landau equation, 

^ = (1 + - m., mM^,t) + K^-^. (2) 

The CGLE has been extensively investigated in the past for its rich equilibrium states and for 
their applications to a variety of physical situations 0,14^. In the presence of time delay it 
is not meaningful to take the continuum limit since a would also make the propagation 
time tend to zero. One can instead take the so called thermodynamic limit by letting N ^ oo 
and the chain length L — oo, while keeping a = L/N fixed. In such a limit the time 
delay parameter remains finite and meaningful. The thermodynamic limit also preserves 
the discrete (in space) character of the evolution equations and therefore the corresponding 
equilibrium states can differ in character from those obtained from a continuum limit. We 
will discuss some of these differences in the next section where we derive and analyze the 
dispersion relation for plane waves for our discrete set of equations (P) . 

III. DISPERSION RELATION FOR PLANE WAVES 

We seek plane wave solutions of Eq. (jT)) of the form, 

where a is the distance between any two adjacent oscillators and k is the wave number such 
that —7c < ka < TT. The values of ka, which define the phase difference between adjoining 
oscillators, are further constrained by the periodicity conditions inherent in a closed chain 
configuration. Since ipN+i must be identical to ipi and i/jq must be identical to ipN, we must 
satisfy the condition e*'^*'" = 1. This implies that Nka = 2mn, m = 0,1, . . . , N — 1, that is, 

ka = m—, m = 0,1, . . . ,N - 1. (4) 
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These discrete values of ka are one of the defining properties of the various phaselocked 
states of the coupled set of oscillators. Further characteristics of these states are described 
by a dispersion relation which can be obtained by substituting Q in Eq. ((T)) to give, iu = 
l + iu;Q-R'^ + 2K [cos{ka)e-'^^ - 1] . This leads to 

u = uq — 2K sm{u!T) cos{ka), (5) 
= 1-2K + 2K cos(cjr) cos(A;a). (6) 

Equations (0) and (jH)) which constitute the dispersion relation for plane waves of our model 
equations can be compared to the corresponding dispersion relation obtained from the 
CGLE, where we can substitute ip{x,t) = RcGL^xp^ikx — iucGit) in © to get, 

l^CGL = ^0, (7) 
RloL = 1 - K{ka? (8) 

where we have replaced K by its scaled value. Equations (|7j) and (jH)) can also be recovered 
from (0) and (jH)) by setting r = and taking the long wavelength limit of fca <^ 1. As we 
can see there are significant differences in the two dispersion relations. While (0) and © are 
valid for any arbitrary value of A^, ((Tj) and (jH} are strictly valid only in the limit of — > oo. 
The spectrum of ka values are therefore continuous for the latter whereas in the former they 
assume discrete values that are determined by the magnitude of A^. For r = 0, Eqs. © and 
((Tj) become identical, but Eq. (0) reduces to, 

R'^ = 1-2K + 2Kcos{ka) = 1 - AK sm'^{ka/2). (9) 

Using i?^ > as the defining condition for the possibility of having a phaselocked state, we 
see that this region is considerably reduced in the case of the CGLE as compared to the 
discrete model equations for any given value of the coupling parameter K. For example, at 
K = 1/4, the discrete model allows all modes from to vr to exist, whereas the continuum 
model has an upper cutoff at ka = 2.0. For K > 1/4, cutoff regions also appear for the 
discrete model and are given by the expression, 

= cos"^(l - 1/2K) <ka< f2 = 271- cos"^(l - 1/2K) 

As can be seen from the above expression, for K > 1/4 the anti-phase locked state {ka = vr) 
is now no longer a permitted state. With the introduction of time delay, the existence 
region is not only a function of ka (for a given value of K) but also depends on uj which 
satisfies the transcendental equation This brings about a significant modification of the 
existence domains and raises the interesting possibility of enabling many of the previously 
forbidden modes to exist. Since the dispersion relation now has a transcendental character 
it also introduces additional branches of collective oscillations. To illustrate some of these 
modifications we consider a few simple examples. Let C = cos{ujt). By choosing a value for 
C and in turn defining S = sin(ti;r), we can immediately get an expression for uj, and the 
resultant equation for R^ defines the stability regions. 

Case (i): First let cos(a;r) = 1. This gives uj = luq, and the square of the amplitude as 
Ai = 1 — 2K + 2K cos{ka) . This also leads us to r = 2mT/uQ, n = 0,1,2, .. .. This case also 
includes the special case of r = 0. Since Ai > for the existence of the solutions, it is easily 
seen that the bounding region in [K, ka) plane is defined by 

fi<ka< f2 = 2n- /i 
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FIG. 1: Allowed (unshaded) and forbidden (shaded) wave modes in the presence of time delay for 
various values of cos(u;r) 



where /i = cos~^(l — 1/2K) is the forbidden region. Since /i does not intersect with the 
K axis at any finite K, the region exists for all -ft" > 1/4. In particular the in-phase state 
exists for all A" > 0, and the anti-phase locked state {ka = vr) no longer exists if > 1/4. 
This forbidden region is illustrated in Fig. ^a). 

Case (ii): Let cos(c<jr) = 1/2. This gives the dispersion relation uj = uoq — v^cos(/ca). 
The corresponding delays are defined by r = (2n±l/3)7r/[u;o — -\/3-ft'cos(/ca)], n = 0, 1, 2, ... . 
The square of the amplitude is given by = 1 — 2K + K cos( fca) . Let gi = cos~^ (2 — 1/K). 
The curve gi assumes values for K between 1/3 and 1. If ii' > 1, ^2 < for any ka. So all 
the modes are forbidden if K > 1, and the region between K = 1/3 and iiT = 1 as defined 
by 

gi<ka< g2 = 2n - gi. 

is also forbidden. In particular the in-phase solutions do not exist if A' > 1 and the anti-phase 
solutions do not exist if K > 1/3. This region is illustrated in Fig. ^b). 

Case (in): Let cos(u;r) = —1/2. This gives the dispersion relation u = uiq + v^cos(/ca). 
The corresponding delays are defined by r = (2n±l/3)7r/[ci;o + \/3A'cos(fca)], n = 0, 1, 2, ... . 
The square of the amplitude is given by A3 = 1 — 2K — K cos{ka) . Let hi = cos~^ {1 / K — 2) . 
If > 1, ^3 < for all ka. So > 1 is the forbidden region. Between K = 1/3 and K = 1 
the region defined by 

hi > ka > h2 = 2-n — hi 
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FIG. 2: Dispersion relation between allowed wave numbers and the corresponding frequency shown 
as r is gradually increased. K = 1 and loq = 10. A range of r values is forbidden. 

is forbidden. The inequality signs are reversed in this case because the curvature of hi is 
different from that of gi or /i. In particular the in-phase solutions are forbidden if i^' > 1/3 
and the anti-phase solutions are forbidden if > 1. This region is illustrated in Fig. E^c). 

Case (iv): Let cos(u;r) = —1. This results in the dispersion relation oj = ujq just as in 
case (i). This case corresponds to r = (2n7r + 1)11/^0,71 = 0,1,2,.... The square of the 
amplitude is given by = 1 — 2K — 2K cos{ka). Since = 1 at ka = vr, in contrast to the 
three previous cases, the anti-phase locked solutions exist for all values of K. At ka = 0, 
A4 = 1 — 4K. So the in-phase solutions are forbidden for all iiT > 1/4. The forbidden regions 
when K > 1/4 are defined by 

ii > ka > i2 = 2-71 — ii, 

where ii = cos~^(l/2i^ ~ !)• This region is illustrated in Fig. ^d). In fact we can derive 
a general expression for the existence curves by choosing any arbitrary value for C. This 
in turn defines the frequencies a.s u = Uq — 2KS cos{ka). After some simple algebra, the 
corresponding delays are given as 

T = {2mr + cos-^C)/[tuo-2KScos{ka)], if S>0, 

T = {2mT - cos"i C)/[ujo - 2KScos{ka)], if S < 0. 

The curve that defines the boundary of the forbidden region is given by ji = cos~^(l/5' — 
1/2KS) when S" 7^ 0. The case of 5 = is the same as case (i) studied above. 

Apart from these special cases, the general existence regions are complicated functions 
of ka, K and r. They need to be determined numerically by simultaneous solution of 
equations (0) and (0). To demonstrate the constraint imposed by (0) we have plotted its 
solution {uo vs. ka) for various values of r and for a fixed value of c^o and K in Fig. El^a). 
When r = 0, the allowed range of r is given by abs(/ca)< cos~^(l — l/2i^). So for K = 1, 
the phase locked patterns that have wave numbers less than vr/S are allowed, and all of 
them have an identical frequency. As r is increased the frequency of oscillation decreases 
for small r, and the dispersion relation acquires a nonlinear parabolic character. As r is 
further increased, depending on the actual value of K, there are bands in r values where 
no modes exist. The shrinking and disappearance of the dispersion curve at fca = beyond 
r = 0.125 up to r = 0.2 in Fig. |2fa)-(b) illustrate this phenomenon. One also notices from 
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Fig. |2tb) that at higher values of r the dispersion curves become discontinuous and have 
bands of forbidden ka regions. 



IV. STABILITY OF PHASELOCKED SOLUTIONS 



We now find the stabihty of the equihbrium phaselocked solutions discussed in the pre- 
vious section by carrying out a linear perturbation analysis. Let, 



J(jka) 



(10) 



where = 0, 1, ■ ■ ■ , — 1. Substitution of (fTUj) in Eq.(IT)), yields in the lowest order the 
dispersion relation discussed in the previous section. In the next order, where we retain 
terms that are linear in the perturbation amplitude, we observe that 



I ^fc I ^fc = h 

Using the above, we obtain the equation 
duj{t) 



2Iiiu, {t) + Rle^^--% it) 



^i{jka) _|_ 



(11) 



dt 



(1 + iujo - 2Rl - 2K) Uj it) - Rle^''^^% (t) 



(12) 



+K 



u,+i {t - r) e^^^" + (t-r). 



-ika 



and, taking its complex conjugate, 

duj{t) 



dt 



l-iujQ- 2Rl - 2K) Uj (t) - Rle-^'^^^^Uj (t) 



(13) 



+K Uj+i {t - r) e-''"' + Uj^i {t - r) 



Multiply Eqs. (fT^ and (fT^ term- by-term by e'^^^'^"'\ make use of the identities 

Uj±i (t - t) e^^^'e'^^'^''^ = Uj±i (t - r) 

and 



J(j±l)qa^±i(k-q)a 



Uj±i {t - t) e^'^'^e'^^''^^ = Uj±i {t - t) e 



(14) 



(15) 



and sum over j = 0, 1, 2, ■ ■ ■ , — 1. Introducing adjoint amplitudes Wg (t) and Wg (t) by the 
definitions 



N-l 

E 

j=0 



[Wg (t) , Wg (t)] = (^) > (^)] e 



i{jqa) 



(16) 



we obtain the set of coupled equations 
dwq it) 



dt 



(1 + lujQ - 2Rl - 2K) Wg it) - Rle^^'^^'wg (t) 
+2K cos [(A; — q) a] Wg {t — r) 



(17) 
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and 

= (1 - ^coo - 2Rl - 2K) w, (t) - Rle-'^^^'w, (t) 



(18) 



+2K cos [{k + q) a] Wg (t - r) . 
In order to perform the stability analysis, we assume solutions of the form 

K (t) , Wg (t)] = [ce'^"', ce-'^^'] (19) 

which yield the set of coupled equations 

C+c + Rlc = 0, (20) 
Rlc + C„c = 0. (21) 

In these equations, the quantities C± are given by 



C+ = X 



l±i{iUo- LJk) - 2Ri - 2K\ - 2K cos [(A; t q) a] e^'"^"^ e'^^ . (22) 

The eigenvalue equation is obtained from the determinantal condition for Eqs. ()20|) and ()2ip. 
namely, 

-Rl = Q. (23) 

This can be expanded to be written in the form of the following characteristic equation, 

+ (ai + a2)A + (aiOa - i?^) = 0, (24) 

where ai = 2R^ - I + 2K - i{uo - u) - 2K cos[{k - q)a]e-^^+'^^^ , a2 = 2R^ - 1 + 2K + 
i{uJo -u)~2K cos[{k + q)a]e~^^~'^^\ 

It should be noted that the perturbation wave numbers q are once again a discrete set 
and from the periodicity requirement they obey the relation, 
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qa = m — , m = 0, 1, . . . , — 1. 

Thus in our stability analysis any pattern corresponding to a given value of ka, we need to 
examine the eigenvalues of Eq. (j24|) at each of the above permitted values of qa. We now 
proceed to discuss the stability of the various phaselocked patterns both in the absence and 
presence of time delay. 

A. Stability of phaselocked patterns in the absence of delay 

In the absence of time delay, the eigenvalue equation (|24j) can be solved analytically to 
give. 



X = 1-2R'^ -2K + 2K cos{ka) cos(ga) ± ^^47^2 sm^{ka) sm^{qa) + (25) 

The real parts of the eigenvalues of Eq. ()24|) will be negative if ai + 02 > and aia2 — R^ > 
simultaneously. The first of these conditions can be simplified to give 

ai + a2 = 2[l-2/s:{l-cos(A;a)[2-cos(ga)]}] > 0, (26) 

and the second condition can be simplified to give 

0x02 — -R"^ = [1 — cos(ga)]4_ft' K ^Acos^{ka) — 2cos{ka) — [1 + cos(ga)]| + cos(A;a) > 0. 

(27) 

In the following, we use these two conditions, or, in the simplest cases, the eigenvalue 
equation itself to determine the stability. 
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1. In-phase patterns (k = mode) 



From the previous section the k = solution is given by cu = luq and R = 1. Substituting 
these in Eq. above, and setting A; = 0, we get, 



A 



-2K[1 - cos(ga)] < 0, 
-2 - 2K[1 - cos(ga)] < 0, 



Thus the plane wave solution "^jit) = e*'^°*, which is nothing but an in-phase locked solution 
of the coupled identical oscillators, is stable for r = for all values of Uq. 



2. Anti-phase patterns (ka = tt mode) 

The equilibrium pattern in this case consists of adjacent oscillators remaining vr out- 
of-phase at all times and oscillating with the same frequency uq- The amplitude of the 
oscillations is given by, 

= 1 _ 4K. 

As discussed in the previous section, the condition > permits the existence of this 
mode in the region of i^' < 1/4. However as we will see below, using conditions (j26j) and 
(j27p . these permitted anti-phase states are linearly unstable for any arbitrary number of 
oscillators. Inserting ka = tt in (j26|) and (j27p the conditions simplify to, 

K < min- ^ , (28) 

2[3-cos(ga)]' ^ ^ 

and 

K > max -— . (29) 

5 — cos(ga) 

where the qa values are to be determined as per the prescription discussed in the previous 
section. As simple examples, for N = 2 the only permitted values of qa are and vr and 
from ()28p and ()29|) above, we arrive at the requirement that K < 1/8 and K > 1/4 which 
is not possible. Hence the mode is unstable. Such conditions were previously obtained by 
Aronson et al. [16] in their investigation of the collective states of two coupled limit cycle 
oscillators. In fact, the argument can be extended to any arbitrary even N since 

5 — cos(ga) < 2[3 — cos(ga)] 

for all values of qa and hence the conditions (|28|) and (|29|) cannot be satisfied simultaneously 
for any value of K. This implies that the anti-phase states, characterized by ka = vr, are 
unstable for any arbitrary value of even N. 



3. Other phaselocked patterns (k > 0) 

Another interesting phaselocked pattern is the k = tt/2 mode. It can be shown by using 
the conditions discussed above that this pattern is also always unstable. However, a large 
number of other modes whose wave numbers are close to either or 2tt are likely to be 
stable. They coexist with the in-phase stable solutions {ka = 0). The relations and 
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FIG. 3: The contours of ai + 02 = (thin hnes), and that of aia2 — i?^ = (thick hnes). The thin 
+ and — signs indicate the sign of ai + 02 in the connected regions bounded by the thin curves, 
and the thick signs indicate the sign of 0102 — R'^ in those bounded by the thick curves. 




FIG. 4: The stabihty region of the phaselocked states (shaded region which is below K = K*) 
plotted in (K, ka) space. The left half of the dashed curve is /i and the right half is /2 as drawn 
in Fig. UJa). 

(12 7|) provide sufficient conditions to find the stability of any given mode. The contours of 
ai + 02 = 0, and 0102 — -R^ = 0, are given, respectively, by 

qa = {cos"^[2 - {2K - l)/2KC],2n - cos"^[2 - {2K - 1)/2KC]] , 
qa = {cos~^ [C{l/K - 2 + 4C) - 1] , 27r - cos"^ [C{l/K - 2 + 4C) - 1]} , 
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FIG. 5: Multiple stable phaselocked patterns for N = f>Q aX loq = 10 and K = 0.4 in the absence 
of time delay. The real parts of the complex vector Zj{t) are plotted in gray scale. 



where C = cos(fca). These curves are plotted in Fig. |21for K = 0.4. We see that there is a 
large range of ka values that satisfies the conditions and (j2Z|)- For a given number 
of oscillators if the ka values fall in this range then the corresponding phaselocked patterns 
are stable. The x-axis intersections of the curve 0102 — = that are closer to and 27r 
provide boundaries of K below which the modes are stable. By setting ga = in (P7j) . we 
obtain the condition on K for stability as a function of ka: 

— cos{ka) _ , . 

^ 2[2cos2(A;a) -cos(A;a) - 1] " ' ^ ^ 

This stability region is plotted in Fig. |3] In the limit of — 00, a continuous range of 
ka modes are accessible, and the system truly possesses infinitely many stable phaselocked 
states when K < K*. The stability of these phaselocked solutions is one important result 
of our paper. Note that all phaselocked patterns with wave numbers between tt/2 and 
37r/2 are unstable. We also note from Fig. |3]that for a given value of K it is possible to 
have more than one stable state corresponding to different values of ka that lie in the two 
stable regions. In the limit oi N ^ 00 such a multistability phenomenon can occur over 
a continuous range of ka values spanning the stable regions. A numerical example of the 
multistabihty of some of these modes is illustrated in Fig. El for A^ = 50 oscillators at a fixed 
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FIG. 6: The in-phase {ka = 0) locked frequencies and amplitudes represented by Eq. ^ at K = 10, 
and ujQ = 10. The dashed portions represent the unstable branches for the case of = 2 and the 
continuous lines represent the stable branches. 

value of K = 0.4. By giving initial conditions close to the modes ka = 27r/50, 47r/50, and 
677/50, the corresponding phaselocked solutions are realized. 

Now we address the following question: What is the minimum number of oscillators for 
which a second phaselocked pattern (the first one being in-phase) can become stable? Since 
the permissible ka values are related to the number of oscillators N in an inverse fashion due 
to the periodicity requirement (see Eq. (@))) we need to determine the maximum stable ka 
in order to find the minimum critical N. Since K = K* makes an intersection at ka = tt/2 
and ka = 3tt/2, the second stable phaselocked pattern that could emerge at this point will 
have 



Hence > 4 for a second phaselocked pattern to ever become stable. Thus for N = 2 
or N = 3 coupled identical oscillators we will only have in-phase oscillations as stable 
oscillations and there will not be any other stable phaselocked states. The minimum number 
of oscillators necessary for multistable behavior is four. Thus = 4 is a critical number. 

B. Phaselocked patterns for finite time delay 

Now we study the synchronized patterns that could become stabilized by the presence of 
time delay. We noted before that, in the absence of time delay, the in-phase oscillations are 
always stable, and the anti-phase oscillations are always unstable for any positive finite K. 
Time delay changes the scenario and could make each of these branches stable in certain 
ranges of r. In addition, the number of multiple branches of each of these two oscillations 
increases with r. Equations (jSlEI) define the amplitudes and frequencies of these synchronized 
states, and the characteristic equation (j^^ determines their linear stability. 

We begin with a simple example of two coupled oscillators, N = 2, which can have just two 
phaselocked states, namely in-phase states with ka = and anti-phase states with ka = vr. 
We evaluate the eigenvalues for these two cases for different values of r by numerically solving 
Eq. ()24|) for the permitted values of qa and also confirm the existence of the modes by an 
actual integration of the original equations. We show the stable (continuous) and unstable 
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FIG. 7: The anti-phase {ka = vr) locked frequencies and ampUtudes represented by Eq. at 
K = 10, and loq = 10. The dashed portions represent the unstable branches for the case of = 2 
and the continuous lines represent the stable branches. 
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FIG. 8: Real part of the eigenvalues as a function of qa for = 10 at cjo = 10, K = 10, ka = tt, 
and T = 0.2. 

(dashed lines) branches of these two solutions in Figs. IHl and |7| A stable in-phase branch 
emerges from r = and its amplitude becomes while its frequency is finite showing a 
supercritical Hopf bifurcation. The first branch of the anti-phase locked state emerges from 
zero again in a supercritical Hopf. As r is increased multiple Hopf bifurcation points are 
seen. Such stability of the in-phase and anti-phase oscillations were earlier studied by us 
experimentally using coupled electronic oscillators j4o|. 

As a second illustration, we consider the case of = 10, r = tt/cjo (i.e. cos(ci;r) = —1.0), 
whose equilibrium existence domain we have discussed earlier as case (iv) in Section UTTl In 
particular we have noted that the ka = n state is a permitted equilibrium state. To examine 
its stability we observe that we need to examine the eigenvalues for qa = m^, m = 
0, 1, . . . , — 1. In Fig. ISl we have plotted the real part of the eigenvalues for all these 
perturbation wave numbers. We see that they are all negative, which indicates stability of 
the pattern. In fact, in this case, due to the symmetry of the characteristic equation we can 
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FIG. 9: (a) Critical boundaries below which stable phaselocked patterns exist are shown for the 
case of cos(a;T) = 1 as ^ oo. (b) Stability region of phaselocked patterns for cos(a;r) = —1 as 
iV ^ oo. 
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FIG. 10: Anti-phase locked state of identical oscillators for = 50 at X = 10 and r = 0.3. 



predict stability for all higher values of as well with the eigenvalues of all additional qa 
values falling on the dotted line shown in the figure. The existence regions for this case were 
shown in Fig. ^d). By considering this case for all the possible values of qa between and 
2tt, we are essentially determining the stability of this case in the infinite oscillator limit. We 
show our numerical results for this case in Fig. El It is obtained by tracking the eigenvalue 
transitions for qa between and 27r. This phase diagram is strictly true for N oo. As can 
be seen, finite r increases the number of stable modes for any given value of K. We finally 
show in Fig. a numerical example of the anti-phase oscillations using N = 50 oscillators. 



15 



V. AMPLITUDE DEATH FOR FINITE TIME DELAY 



As is well known, a system of identical oscillators that are globally coupled do not have 
an amplitude death state jl^ but the presence of time delay in the coupling can bring about 
such a collective state. This phenomenon was first pointed out for two coupled oscillators in 
[35I I and also generalized to globally coupled oscillators js^. A similar statement also holds 
for locally coupled oscillators as will be shown. The death state in the presence of time delay 
had earlier been numerically confirmed for a finite number of locally coupled oscillators |,35| . 
However no systematic study of the dependence of the death island regions on the magnitude 
of the time delay and the number of oscillators has so far been carried out for systems of 
locally coupled oscillators. In this section we address this issue and study the stability of 
the origin for our discrete model equation ((T)). Amplitude death state is characterized by 
R = 0. To derive an appropriate set of eigenvalue equations for determining the stability 
of this state, we carry out a linear perturbation analysis about the origin. Substituting 
{ipj = + ipj) in Eq. (Q), and discarding the nonlinearities, we get, 

^ = (1 + iuJo)i', + - r) - i^M + - r) - ^,{t)], (31) 

where j = 1, . . . , A^ with periodic boundary conditions. By letting iplt) oc e*'^*, the eigenvalue 
matrix in circulant form is obtained. The determinant of this matrix is written as 

N 

Y[(^X + 2K -l-icuo- Ke-^^Uj - Ke-^^Uf'^) = 0, 
i=i 

where Uj = e^2^(J-i)/^ are the A^*^ roots of unity. But Uj + Uf'^ = Uj + Uj^ = 2cos[(j - 
l)27r/A^]. So the above equation takes the form of 

N 

l[(^X + 2K -l-iuJo-2K cos[(j - l)27r/A^]e"^") = 0. (32) 

The complete set of eigenvalue equations includes the second set obtained by considering 
the conjugate equation of the above. Note that for r = the above eigenvalue equation 
fj32|l always admits at least one unstable eigenvalue, namely A = 1 + iujQ. Hence identical 
oscillators that are locally coupled cannot have an amplitude death state in the absence 
of time delay. We will now determine the amplitude death regions for finite values of r. 
We define a factor Rj = 2cos[(j — l)27r/A^] that we will use in the critical curves derived 
below. If the number of oscillators is a multiple of 4, there are some eigenvalue equations 
that emerge without a dependence on r, when Rj = 0. For example, consider the case of 
A^ = 4 and j = 2,4. Then, the eigenvalue equation becomes: A = 1 — 2K ± iuq. For this 
equation, the only criticality is given by A' = 1/2. The stable region lies on the side of the 
parameter space that obeys K > 1/2. For other values of Rj, the death island boundaries 
can be derived by setting the real part of the eigenvalue to zero, and appropriately choosing 
the signs of the multiple curves that result. The analysis is similar to the one we presented 
in the treatment of globally coupled oscillators |3^ i3^] , and here we simply provide the final 
expressions for the critical curves in (r, K) plane: 

2mr-cos-'^[{2K~l)/KRj] D ^ f) 



WO' 
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{ 27nn+cos-'^[{2K-l)/KRj] /? \ f) 
(2m+l)7r+cos-M(2i^-l)/A-|fl, |] D ^ n ^ ^ 

n and m are whole numbers. We now determine some useful bounds on K for ordering 
and finding the degeneracies of the critical curves. The argument of inverse cosine functions 
and the square root term in the denominators of the above expressions impose the following 
bounds on K: 

K > 1/4, for Rj=i, (35) 
— - — <K < — - — , for i?,>i. (36) 

The sign of the derivative with respect to r of the real part of the eigenvalue is determined 
by the term — Im(A)ii'i?j sin[Im(A)r], and after some algebra it is expressed as 

for any uq, 

dRe{X) 



dr 




if uo > ^Km] - {2K - 1)2, (37) 

Re(A)=0 > 0, on Tfe, if UJo< K^] - {2K - 1)^ . 

The first condition involving uo in the above relation leads to the following bounds on K: 

< (l + cu2)/4, for (38) 



2 + i^V^+^ 2-i?^/l+^ 

4^^2 ' ^ < ^3^2 ' for ^i>i- (39) 

The ranges of K imposed by the relations and (IHUj) are mutually exclusive. Hence, on 
the curve Tb the only boundary across which an eigenvalue pair makes a transition to the 
negative eigenvalue plane occurs when Rj=i, and across the curves occurring when Rj^i, the 
eigenvalue pair makes a transition to the positive plane. So rf,(0, K) at Rj=i always remains 
as the left hand side boundary of the death island. In order to see the degeneracy among 
the curves, note that 

Rk = Rj, if k = N + 2-j. 

Also the sign of Rj does not play a role in distinguishing the critical curves. These two 
properties are responsible for a reduction of the number of the actual distinct boundaries 
of the death islands. These degeneracies can be framed into two cases. First, when is 
odd. The total number of distinct values of Rj is (A^ + l)/2. -Ri(= 2) is the maximum 
of all the Rj, and is non-degenerate. And Ri, where i = 2, . . . , (A^ + l)/2 are the other 
distinct values, whose values are identically equal to Rk, where k = N + 2 — i. Among 
the latter Rj values, the maximum negative value occurs at j = {N + l)/2. The curves 
corresponding to Rj > and Rj < form the boundaries for two different death islands, as 
can be seen from the indexes in Eqs. ()33H34|) . A further ordering of the curves must be done 
by a numerical plotting of the curves. The ordering reveals that for A^ < 13, the first death 
island is bounded by the curves 
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Th{0,K) at Rj=i, 

n{0,K) at Rj=(^N+i)/2, and 

Ta{0,K) at Rj=(N+l)/2, 

Ta{l,K) at Rj=i. 
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FIG. 11: Death islands at ujq = 10. All even number of oscillators have a single death island region 
that is independent of the number of oscillators. The odd number of oscillators are bounded by 
four curves when N < 13, and two curves otherwise. These two curves merge in the infinite limit 
with the curves that represent the even number of oscillators. 



The curves 71 and 62 form, respectively, the left and the right boundaries of the death island 
region. And and 62 form the bottom two curves. The existence range along K of Si, across 
which the eigenvalues make a transition to the left half plane, increases with increasing A^. 
In fact it intersects with 71. So this provides the first boundary across which stability is 
lost. This occurs for > 15, when the death island is bounded by the two curves: 

7i : Th{0,K) at Rj=i, and 

61 : Ta{0,K) at Rj={N+l)/2- 

Second, when A^ is even, many more curves become degenerate. The number of distinct 
value of Rj is N/2 + 1. Among these distinct values, the magnitudes of pairs of them can 
become identical. If A^ is divisible by 4, such pairs are N/4 + 1 in number, and (A^ + 2)/4 
otherwise. Every positive Rj has its negative counterpart. The maximum value of Rj = 2 
that occurs when j = 1, and the next maximum value is at Rj=2- Hence the equations 
()33II34|) can be simplified to 

Te n, K) = , " Rj=i J, (40) 

, m7r + cos-i[(2ir-l)/A^i?,] ^ 

LUo - ^Km] - {2K - 1)2 

where J = N/4 + 1 if A^ is divisible by 4, and (A^ + 2)/4 otherwise. The actual ordering, 
however, reveals that the death island boundaries are given by Tc{0, K) at Rj=i and rrf(l, K) 
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FIG. 12: Amplitude death of identical oscillators for = 50 at X = 10 and r = 0.1. 



at Rj=i, which are identical, respectively, to 6i and 71. These curves are plotted in Fig. ITTl 
As is seen for even there is a single death region. In the case of A^ odd the boundary of 
the death region depends on the value of A^. As A^ increases, the area of the death region 
decreases. As A^ ^ 00 the area of the death island for A^ odd decreases and approaches, 
as a limit, the boundary for A^ even. For all A^, the intersections of 71 and 61 or those 
of 72 and 62 occur for i^' > 1/2. So the delay-independent eigenvalue equations that are 
mentioned earlier do not contribute to the death island boundaries. The differences in the 
death island boundaries for even and odd numbered oscillators can be traced primarily to 
the behavior of the eigenvalues of the lowest permitted perturbation wave numbers. For an 
even number A^ of oscillators the smallest perturbation mode is qa = vr. The values of the 
real parts of the eigenvalues corresponding to this mode are close in their magnitude to those 
corresponding to the qa = perturbation mode. Across the right hand side boundary of 
the death island region, the qa = ir mode grows positive and the system emerges out of the 
death region with an anti-phase state. When A^ is odd, however, the smallest perturbation 
mode is qa = IT — 7c/N which is more heavily damped than the qa = mode. So the death 
region continues to exist for larger r values. Ultimately, the second eigenvalue branch of 
the qa = mode (which exists due to the transcendental nature of the eigenvalue equation) 
grows and the system emerges across the boundary with an in-phase state of a different 
frequency. As A^ becomes large, the smallest perturbation mode for the A^ odd case gets 
closer to n and the death island boundaries of the two cases, as seen in Fig. become 
indistinguishable. We have independently verified the death regions depicted in Fig. [TTl 
including their interesting dependence on A^, by a direct numerical solution of (0) over the 
specified range of parameters. We also illustrate in Fig. an example of the amplitude 
death state for a chain of A^ = 50 oscillators. 
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VI. DISCUSSION AND CONCLUSIONS 



We have studied the existence and stabihty of phaselocked patterns and amphtude death 
states in a closed chain of delay coupled identical limit cycle oscillators that are near a 
supercritical Hopf bifurcation. The coupling is limited to nearest neighbors and is linear. 
The coupled oscillators are modeled by a set of discrete dynamical equations. Using the 
method of plane waves we have analyzed these equations and obtained a general dispersion 
relation to delineate the existence regions of equilibrium phaselocked patterns. We have also 
studied the stability of these states by carrying out a linear perturbation analysis around 
their equilibria. Our principal results are in the form of analytic expressions that are valid 
for an arbitrary number of oscillators (including the N ^ oo thermodynamic limit) and 
that can be used in a convenient fashion to identify and or obtain stable equilibrium states 
for a given set of parametric values of time delay, coupling constant, wave number and wave 
frequency. We have carried out such an exercise for a number of illustrative cases both 
with and without the presence of time delay. In the absence of time delay, our analysis 
reveals a number of new phase-locked states close to the in-phase stable state which can 
exist simultaneously with the in-phase state. The minimum number of oscillators for which 
this multirhythmic phenomena can occur is N — 4. Time delay introduces interesting new 
features in the equilibrium and stability scenario. In general we have found that time delay 
expands the range of possible phaselocked patterns and also extends the stability region 
relative to the case of no time delay. The dispersion curves for varying values of time delay 
also display some novel features such as forbidden regions and jumps in the range of allowed 
wave numbers as well as forbidden bands in the space of time delay. 

We have also carried out a detailed analytic and numerical investigation of the existence 
of stable amplitude death states in the closed chain of delay coupled identical limit cycle 
oscillators. The results not only confirm our earlier numerical demonstration of the existence 
of such death states but go beyond to provide a comprehensive picture of the existence 
regions in the parameter space of time delay and coupling strength. The analytic results 
also establish that death island regions exist for any number of oscillators for appropriate 
vales of K, ui and r. In this sense our work provides a generalization of the earlier amplitude 
death related results, that were obtained for globally coupled oscillators to the case of locally 
coupled oscillators. An interesting new result, arising from the local coupling configuration, 
is that the size of a death island is independent of A^ when A^ is even but is a decreasing 
function of A^ when A^ is odd. In other words the death island results for the N = 2 island 
hold good for any arbitrary even number of locally coupled oscillators and constitute the 
minimum size of the death island in the K — t parameter space. This can have interesting 
practical implications. For example in coupled magnetron or laser applications if one is 
seeking to minimize the parametric region where death may occur (and thereby greatly 
diminish the total power output of the system) it is best to select a configuration with an 
even number of devices. At a more fundamental level this "invariance" property which is 
strongly dependent on the symmetry of the system may also have interesting dynamical 
consequences, e.g. in the manner of the collective relaxation toward the death state from 
arbitrary initial conditions. In fact we see some evidence of this differing dynamics in our 
numerical investigations of the time evolution of the system toward the death state. For an 
even number of oscillators we observe a rapid clustering of the oscillators into two distinct 
groups that are tt out of phase. These two giant clusters then slowly pull each other off 
their orbits and spiral toward the origin. For an odd number of oscillators the lack of 
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symmetry appears to prevent this grouping, the phase distribution of the oscillators has 
a greater spread in its distribution and the relaxation dynamics is distinctly different. As 
the number of oscillators increases and the asymmetry gets reduced the difference in the 
dynamical behavior becomes less distinct. In the hmit of — > oo, the difference vanishes 
and the size and shape of the odd island asymptotes to the A^-even island. A fundamental 
understanding of this dynamical behavior and its relation to the symmetry dependence 
emerging from the stability analysis could be an interesting area of future exploration. 
Our results could also be useful in applications where locally coupled configurations are 
employed such as in coupled magnetron devices, coupled laser systems, and neural networks 
as a roadmap for accessing their various collective states. 
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